﻿using System;
using MathNet.Numerics.Distributions;
using MathNet.Numerics.LinearAlgebra.Double;
using MathNet.Numerics.LinearAlgebra.Generic;

namespace Marvin.Optimization
{
    /// <summary>
    /// Implementation of gradient descent.
    /// </summary>
    public class GradientDescent: IOptimizationAlgorithm
    {
        private double _initializationRange;

        #region Implementation of IOptimizationAlgorithm

        public GradientDescent(IGradientFunction gradientFunction, double learningRate = 0.001,
                               double precision = 0.0000001, double initializationRange = 100.0)
        {
            GradientFunction = gradientFunction;
            LearningRate = learningRate;
            Precision = precision;
            InitializationRange = initializationRange; 
        }

        /// <summary>
        /// Gets or sets the learning rate / step size.
        /// </summary>
        /// <value>
        /// The learning rate / step size.
        /// </value>
        public double LearningRate { get; set; }

        /// <summary>
        /// Gets or sets the gradient function.
        /// </summary>
        /// <value>
        /// The gradient function.
        /// </value>
        public IGradientFunction GradientFunction { get; set; }

        /// <summary>
        /// Gets or sets the precision.
        /// </summary>
        /// <value>
        /// The precision.
        /// </value>
        public double Precision { get; set; }
        
        /// <summary>
        /// The initialization vector will be a random position, whereas in each dimension 
        /// the position is between -InitalizationRange and +InitializationRange. 
        /// </summary>
        public double InitializationRange
        {
            get { return _initializationRange; }
            set { _initializationRange = Math.Abs(value); }
        }

        private DenseVector CurrentMinimum { get; set; }

        /// <summary>
        /// Calculates the minimum.
        /// </summary>
        /// <returns>
        /// the position of the minimum in the space defined by the dimensions of the gradient function
        /// </returns>
        public Vector<double> CalculateMinimum()
        {
            InitializeCurrentMinimum(); 
            double delta = double.PositiveInfinity;

            while (delta > Precision)
            {
                var gradientAsArray = GradientFunction.Calculate(CurrentMinimum);
                var gradient = new DenseVector(gradientAsArray);
                var step = LearningRate * gradient; 
                CurrentMinimum = CurrentMinimum - step;
                delta = step.AbsoluteMaximum(); 
            }

            return CurrentMinimum; 
        }

        private void InitializeCurrentMinimum()
        {
            CurrentMinimum = DenseVector.CreateRandom(GradientFunction.GetDimension(),
                                                      new ContinuousUniform(-InitializationRange, InitializationRange));
        }

        #endregion
    }
}
